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ABSTRACT 

We derive the tilt equation governing the inclination of a thin self-gravitating gaseous 
disk subject to low frequency global m = 1 bending perturbations. The disk orbits 
under the influence of a dominant central mass. However, self-gravity can be important 
enough that the disk approaches marginal stability to local axisymmetric perturbations 
(Q ~ 1). The vertical restoring forces due to self-gravity and pressure are evaluated 
correct to first order in the aspect ratio H/r. Thus the effects of bending waves 
are included correct to lowest order in the wave crossing rate (H/r)tt, fi being a 
characteristic disk rotation frequency. 

Both free and forced disturbances are considered. The disk response and precession 
frequency induced by the presence of a binary companion in an orbit with general 
inclination to the unperturbed disk plane are derived. When the degree of warping 
and the inclination are small, it is shown that identical results are obtained if, 
alternatively, perturbation of the disk out of an equilibrium plane coinciding with that 
of the companion is considered, the time averaged potential due to the latter being 
incorporated into the equilibrium potential. The condition for the disk to precess 
approximately like a rigid body with a small degree of warping is found to be that 
the density wave crossing time be significantly shorter than the precession period. 
We consider the effects of the presence of a viscosity which can be characterized with 
the standard a parameterization and find that, to the order we work, the precession 
frequency is unaffected, with any change to the inclination of the slightly warped disk 
occuring at a slower rate. For a <C H/r, effects due to both pressure and self-gravity 
are important in the response, while for a 3> H/r, the response becomes dominated by 
self-gravity with pressure effects becoming negligible. 

As an application of these results, we explore the possibility that the recently observed 
warped disk in the active galaxy NGC 4258 is produced by a binary companion. Our 
results indicate that it may be produced by a companion with a mass comparable 
to or larger than that of the observed disk. If the only source of warping is such a 
companion, a small trailing twist is produced by viscosity for modest a ~ 0.1. The 
outer edge of the disk may also be tidally truncated by the companion. 

Subject headings: accretion, accretion disks — galaxies: kinematics and dynamics - 
galaxies: structure — galaxies: individual (NGC 4258) 
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1. Introduction 

Gaseous disk configurations which are dominated by a central mass and in which internal 
self-gravity is important may occur in different astrophysical contexts such as protostellar disks 
(Larson 1984) or active galactic nuclei (Paczynski 1977). 

The most definitive evidence for a gaseous accretion disk in active galactic nuclei is provided by 
the discovery of megamasers (Claussen, Heiligman & Lo 1984) around the nucleus of the mildly 
active galaxy NGC 4258 (Claussen & Lo 1986). Confirmation of Keplerian rotation was obtained 
from the radio interferometer (Nakai, Inoue & Miyoshi 1993) and VLBI observations (Greenhill 
et al. 1995). Based on the correlation between the spatial locations and radial velocities of the 
masers, Miyoshi et al. (1995) deduced that the masers are located at R ~ 0.13 — 0.26 pc around 
a black hole with a mass 3.6 x 10 7 M® (Watson & Wallin 1994; Maoz 1995). Recently the inner 
part of the disk has been modelled with an advection-dominated flow (Lasota et al. 1996). 

Random deviations from the Keplerian rotation curve at the location of the masers are ~ 3.5 km/ s. 
The lack of systematic deviations provides an upper limit for the disk mass of 4 x 10 6 Mq. The 
disk scale height H < 3 x 10 -4 pc is such that H/R < 2.5 x 10 -3 and the mid-plane temperature 
T c < 10 3 K (Moran et al. 1995). In addition, the high-velocity maser sources have negligible 
acceleration and are not colinear with the low-velocity masers (Greenhill et al. 1995). The radial 
dependence of declination of the red-shifted (with respect to the systemic velocity of NGC 4258) 
high-velocity masers is antisymmetric to that of the blue-shifted high-velocity masers (Miyoshi 
et al. 1995). One scenario to account for these observed properties is that the position angle of 
the rotation axis varies continually with radius by up to 0.2 radians (Herrnstein, Greenhill & 
Moran 1996). This corresponds to introducing a small warp into the disk model. The inclination 
of the local orbital plane with respect to the orbital plane at the innermost region of the disk, g, 
being the ratio of the vertical displacement to the local radius, varies with the disk radius. In at 
least one set of models, small-amplitude variations in the angle between the line of sight and the 
rotation axis are also introduced. This corresponds to introducing a "twist" into the warped disk. 
The magnitude and direction of the twist are not yet well constrained. 

In this paper, we study the dynamics of the warped disk in NGC 4258. The main constraints for 
any theoretical model are: 1) the coherence and the smooth radial dependence in the inclination 
g, 2) the negligible acceleration in the high velocity maser sources, and 3) the nearly perfect 
Keplerian rotation curve. 

We explore the possibility that the observed warp may be produced by a binary companion which 
orbits in a plane inclined to that of the disk. Such non coplanar system could be produced if the 
companion approaches the disk as a result of random gravitational perturbations arising from 
the field stars in the nucleus of NGC 4258 (see § |j). Indeed, in that case there is no a priori 
reason to suppose alignment between its orbital angular momentum vector and that of the matter 
which supplies the disk at large distances. Furthermore, in the evolutionary situation we envisage, 
the disk is not expected to have settled to an equilibrium in which its mid-plane coincides with 



- 4 - 



the companion's orbital plane (see § |1.3j ). Capture of a companion is not expected to result 
systematically in the disruption of the disk. The orbit may have indeed become bound when 
the companion was far away from the disk, with a subsequent decrease of the separation due to 
dynamical friction. 

In the model we shall adopt in this paper, the disk is subject to forces due to internal self-gravity, 
pressure, viscosity, and also to the gravitational forces arising from the central object and the 
perturbing companion on an inclined orbit. Since we shall focus on the secular response of the disk 
(see § [2.2D , the perturbation is the same as that produced by a mass distributed uniformly along 
an inclined ring. For a near Keplerian self-gravitating disk with small viscosity, pressure may in 
principle be as important as self-gravity in determining the response. However, for the model 
we adopt in § |3| for NGC 4258, self-gravity turns out to dominate over pressure and viscosity 
in determining the disk response to the perturbing potential. Before going on to develop the 
tilt equation we use to describe near Keplerian disks, we summarize some general results on the 
dynamics of warped disks. 



1.1. Precession of warped disks 

In a classic paper, Hunter & Toomre (1969) showed that an isolated self-gravitating disk subject 
to a vertical displacement generally precesses differentially, and thus cannot sustain a warped 
configuration. However, differential precession is prevented by gravitational torques from the 
distorted disk itself if the disk has a sharp edge. A clear physical description of how self-gravity 
acts to smooth out differential precession is given by Toomre (1983). A similar process occurs 
if the disk orbits in the external potential due to a companion (Hunter & Toomre 1969) or a 
flattened halo whose equatorial plane is misaligned with the disk plane (Toomre 1983; Dekel & 
Shlosman 1983; Sparke 1984; Sparke & Casertano 1988), provided departure from spherical 
symmetry of the external potential is relatively small. In this case, the disk settles into a discrete 
bending mode (representing a warp) which is referred to as the modified tilt mode because in the 
limit that the external potential is spherically symmetric it reduces to the trivial rigid tilt mode. 
The difference in shape and frequency between the rigid and modified tilt modes is due to the 
fact that, in a non spherically symmetric potential, the disk has to bend to alter the precession 
frequency at each radius so that the rate is everywhere the same (Sparke & Casertano 1988; 
Hofner & Sparke 1994). We note that the potential of an oblate halo can be modeled by that of a 
distant massive fixed ring which is misalinged with the disk plane (Toomre 1983). This situation is 
analogous to that being investigated here in this paper, namely a near Keplerian self-gravitating 
disk dominated by a central mass and a smaller orbiting companion which in a time average sense 
gives the appearance of a massive ring. 

Papaloizou & Terquem (1995) have shown that radial pressure forces in a non self-gravitating 
inviscid accretion disk are also able to smooth out differential precession. This process can 
be effective if the sound crossing time through the disk is much smaller than the precession 
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period. When this condition is satisfied, bending waves are able to propagate through the 
disk sufficiently fast so that the different parts of the disk can "communicate" with each other 
and adjust their precession rate to a constant value. This also happens when viscosity is 
present, but the communication becomes diffusive rather than wave-like when the Shakura & 
Sunyaev (1973) viscosity parameter a significantly exceeds the ratio of disk semi-thickness to 
radius (see Papaloizou & Pringle 1983; Demianski & Ivanov 1997). The theoretical expectation 
of uniform precession of inclined disks under conditions of adequate physical communication has 
been confirmed by the numerical simulations of Larwood et al. (1996). 

When the disk precesses uniformly in an external potential, its precession frequency can be 
calculated from the condition that the disk is stationary in a rotating frame. This condition 
implies that the net torque exerted by the external potential and the Coriolis force is zero (by 
virtue of Newton's third law, the disk does not exert a torque on itself). Kuijken (1991) has used 
this condition to calculate the precession frequency of a self-gravitating disk in a flattened halo 
whose plane of symmetry is misaligned with the disk plane. In the limit of small misalignment, his 
result reduces to that given by perturbation theory (Sparke h Casertano 1988). As expected, his 
results also indicate that when the departure from spherical symmetry of the potential gets too 
large, self-gravity can no longer maintain a uniform precession of the entire disk. In the context 
of simulations of gaseous disks in interacting binary stars that were misaligned with the binary 
orbital plane, Larwood el al. (1996) found that when the tidal potential of the companion is large 
enough, forces due to pressure and viscosity in a non self-gravitating thin disk cannot prevent 
differential precession. 



1.2. Timescale for settling to the warp mode 

We have pointed out above that the disk settles into a (modified tilt) mode which enables it to 
undergo rigid precession. This asymptotic state is possible because bending waves transport away 
the energy associated with the transient response (Toomre 1983). In a purely self-gravitating 
disk, this energy is carried out towards the disk outer edge by bending waves whose group velocity 
approaches zero there. The waves pile up with very short wavelength at the edge where they 
would be expected to dissipate. The timescale for settling to the warp mode is given by the 
characteristic time for density waves to propagate through the disk but the process is slowed down 
in the low surface density regions near the outer edge. Hofner & Sparke (1994) found that when a 
galactic disk is subject to the potential of a misaligned flattened halo, the inner parts settle first 
and the outer parts are unlikely to settle into the warp mode within a Hubble time. 

When pressure is present the waves may be reflected from the edge before they attain arbitrary 
short wavelength there. Some dissipation is then needed for the disk to settle into a state of 
near rigid precession. This can be provided by disk viscosity which is least effective on a global 
disturbance such as near rigid body precession. Other transient disturbances would, as we 
have indicated above, have shorter wavelengths and accordingly would be expected to dissipate 
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significantly more rapidly. In this context note that a warp mode which deviates only slightly 
from a rigid tilt has a wavelength which is much larger than the disk. Larwood et al. (1996) 
considered the disk response due to the perturbation of a companion. In their simulations, of non 
self-gravitating inclined disks with pressure and viscosity, the disk was seen to quickly settle into a 
state of near rigid body precession. The time for this to happen was consistent with about ~ 1/a 
orbital periods, this being the decay time of short wavelength bending waves without self-gravity 
(see Papaloizou & Lin 1994, 1995). 



In this paper, we consider the limiting case of a very thin self-gravitating disk (see § 3T) in which 
bending waves propagate with a speed comparable to that of sound and the settling timescale into 
the warp mode is short compared to the viscous evolution timescale of the entire disk. 



1.3. Timescale for settling to the orbital plane 

Papaloizou & Terquem (1995) have shown that if a non self-gravitating disk perturbed by a 
companion on an inclined orbit evolves towards the symmetry plane of the potential (which does 
not necessarily occur), it does so on a timescale related to the rate at which angular momentum 
is transported to the companion which is expected to be at least the disk viscous timescale. This 
conclusion might be expected since the disk inclination as a whole can change only if angular 
momentum is transferred between different parts of the disk and then between the disk and the 
companion's orbit. As the long term disk evolution is driven through the action of viscous stress, 
adjustment should occur so that angular momentum transfer to the orbit, and hence changes to 
disk inclination occur on that timescale or longer. Such an adjustment can occur for example 
through evolution of the location of the disk outer edge away from the companion's orbit. 

This was seen in the numerical simulations of Larwood et al. (1996) and Larwood (1997), in which 
disk inclination was ultimately seen to evolve much more slowly than the precession rate. We 



also show that this is expected to occur for a mildly warped self-gravitating disk in section [2. 4. 6 
below. Thus we do not necessarily expect the inclination of the disk to have evolved so that it 
coincides with the orbital plane. 

Results of numerical simulations of the dynamical evolution of rings highly inclined to the 
symmetry plane of some galactic potential (Katz & Rix 1992; Christodoulou et al. 1992; 
Christodoulou & Tohline 1993) are also in agreement with the above conclusions. However, we 
emphasize that they only hold if the disk is able to find a state in which differential precession 
can be controlled by internal self-gravity and pressure. Settling toward a preferred orientation is 
much faster if differential precession is so important that it cannot be controlled by other forces 
(Steiman-Cameron & Durisen 1988). This happens, for example, when the external potential 
is such that the precession frequency is always not very small compared to the orbital angular 
frequency. However, this is not the case in the near Keplerian disks we consider. There the 
precession frequency is relatively very small at all inclinations. 
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1.4. Plan of the paper 

The plan of the paper is as follows. The response of the disk to the companion perturbing 
potential is analyzed in § |2[ We first give the basic equations and the expression for the time 
averaged perturbing potential we later use. We have taken the latter to be due to a companion 
on an inclined orbit, but any m = 1 perturbing potential proportional to the coordinate z along 
the disk rotation axis could equally well be considered. We then describe the disk equilibrium 
structure, linearize the basic equations and derive the tilt equation governing global variations in 
the disk inclination taking into account self-gravity and pressure. All effects contributing to low 
frequency bending wave propagation across the disk are included to leading order. 

The disk response and precession frequency are derived using two different approaches, which are 
shown to be equivalent in the limit of small inclination of the companion orbit. We consider i) the 
excitation of small amplitude low-frequency bending modes, considered as perturbations of the 
disk away from the companion orbital plane. Here the time averaged potential of the companion 
is included in the equilibrium state. We also consider ii) the zero frequency response of a disk 
forced by the time averaged potential due to a companion in a general inclined orbit. Here the 
companion is not included in the equilibrium state. We go on to discuss the effect of viscosity, 
showing how this has the effect of reducing the importance of the pressure forces in comparison to 
those due to self-gravity once the viscosity parameter a significantly exceeds H/r. We also show 
that when the disk warping is small, the disk inclination evolves on a long timescale. 

In § |H we apply these results to NGC 4258. We derive the disk parameters from the observations 
and present the numerical results. We find that, if the disk is self-gravitating, the observed 
warp may be produced by a companion of comparable mass to that of the disk, with an orbital 
separation of about one and a half times the observed outer radius of the disk. The presence of the 
companion may also cause the disk truncation. We also show that viscous forces induce a twist 
in the orientation of the spin axes of annuli at different radii. For the estimated disk parameters, 
only a small trailing twist is expected to be associated with the warp if a perturbing companion 
alone produces it. At the present moment, there is no strong constraints on the magnitude of the 
twist. But an accurate measure of it can provide useful information on the non conservative forces 
acting in the disk. 

We summarize and discuss our results in § [|. 

2. Disk response 
2.1. Basic equations 

The dynamics of the disk is described by the equation of motion 
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^ + ( v .V)v = --VP-Vf + f v , (1) 

at p 

and the equation of continuity 

^ + V -(pv) = 0, (2) 

where P is the pressure, p the density, v the flow velocity and ^ the total gravitational potential. 
We allow for the possible presence of a viscous force per unit mass f v but we shall assume that it 
does not affect the undistorted axisymmetric disk so that it operates on the perturbed flow only. 

We write * = ^> ex t + ^g, where ^ e xt is the contribution to potential due to the central black hole 
and the perturber and ^g is the contribution to the potential arising from the disk self-gravity 
given by 

*g — —g I ep", ,3, 

Jv |r - r| 

the integral being taken over the disk volume, with r and r' denoting position vectors and G being 
the gravitational constant. 

For our calculations, we adopt a polytropic equation of state P = Kp l+l l n , K and n being the 
polytropic constant and index respectively. Then is a function of r alone, and the sound speed 
is given by c 2 s = dP/dp. 



2.2. Perturbing potential 

We consider the response of the thin self-gravitating disk to a perturbing gravitational potential 
that has a plane of symmetry which is inclined to that of the disk. We suppose that the disk is 
perturbed by a point mass M p orbiting as a binary companion to the point mass M located at the 
disk centre. In our presentation, the disk outer radius is denoted by R and the binary separation 
by D. 

We consider a non rotating Cartesian coordinate system (x, y, z) centered on the central point 
mass M. The z axis is chosen to be along the rotation axis of the unperturbed disk. We shall also 
use the associated cylindrical polar coordinates (r, ip, z). We take the orbit of the perturbing mass 
to be in a plane which has an inclination angle 5 with respect to the (x, y) plane and the line of 
nodes to coincide with the x axis. We denote the position vector of the perturbing mass by D 
with D = IDI. 



The potential ^> ex t arising from both the central and the orbiting masses is given by 
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GM GM V GALr.D 

®ext = — j— r " 1 ^TT + 



r - D I D 3 

The last indirect term accounts for the acceleration of the origin of the coordinate system. We are 
interested in the warping of the disk which is excited by terms in the potential which are odd in z 
and which have azimuthal mode number m = 1 when a Fourier analysis in ip is carried out. We 
consider only the secular term in the potential, that is to say the zero-frequency term, because 
this is the dominant contribution for producing the large-scale warp structure, in which we are 
interested. This term is obtained by taking a time average of the total perturbing potential, and is 
equivalent to replacing the full potential by the one obtained if the perturbing mass is distributed 
uniformly along the orbit. Non zero-frequency terms would be important in determining the disk 
response if it extended much beyond the inner Lindblad (2:1) resonance with the companion. 
Here we suppose that truncation of the disk, through the action of bending and other waves 
excited by interaction with the companion, occurs such that D = 1.5R (see for example Lin & 
Papaloizou 1993), which more or less excludes the resonances required for effective wave excitation 
from the disk. Non secular bending waves propagating into the disk interior have short wavelength 
and do not give a global response. This is supported by the numerical simulations of Larwood et 
al. (1996) who found that the global precession of truncated disks could be accounted for by the 
secular response. 

The term in the Fourier expansion of the potential which is of the required form is given by 



sin ifi /" 27r /" 27r 

= l^T J J [*«rt(r, <f/, z, ujt) - V ext (r, tpf , -z, wt)] sin tf/dtf/, (4) 

where u is the companion's orbital frequency. The parameters we shall use in the numerical 
calculations require a development to eighth order in D^ 1 of this integral: 



3 GM„ 



ext 



4 A 5 / 4 D 3 



rz sin 25 sin ip 



1 + 



ai 



A 



+ 



02 



A 2 



+ 



a 3 



A 3 



0,4, 



r Y 



, (5) 



with ,4= (1 + r 2 / X> 2 ) 2 and 



ai = 2.1875 (2 - 1.5 sin 2 tf) , 

a 2 = 6.0156 (3 -4.5 sin 2 5 + 1.875 sin 4 6) , 

a 3 = 18.3289 (4 -9 sin 2 5 + 7.5 sin 4 5 - 2.1875 sin 6 5) , 

a 4 = 59.2022 (5 -15 sin 2 5 + 18.75 sin 4 5 - 10.9375 sin 6 5 + 2.4609 sin 8 tf) . 
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Since we consider the thin disk limit, additional terms, smaller by factors containing powers of 
z/r, have been neglected. Throughout this paper we shall use the complex potential 

<*t = i l - m ^rf(r)z siniWe* (6) 

defined such that its real part is equal to the physical potential. The function f(r) is the term in 
brackets in (||). 



2.3. Equilibrium structure of the disk 

In the absence of perturbation, the disk is axisymmetric so that in cylindrical coordinates the 
velocity v = (0,rfi,0) and the hydrostatic equilibrium equations are satisfied in the form 

p dr dr 

i£ = -£. (8) 
p Oz Oz 

In the absence of perturbation or orbiting companion, the disk is under the influence of the central 
point mass M so that ^ ex t = —GM/^/r 2 + z 2 . In the present context, we consider thin disks for 
which the radial scale over which physical parameters vary is very much larger than the vertical 
scale height. In the limit where the Toomre parameter Q is of order unity, the disk's self-gravity, 
which may be evaluated as if the disk had zero thickness, is then more important by a factor of 
order r/H than pressure in ([7]). The angular velocity is given by equation (|^) as 

rp or r or 

1 /2 

where Q K = (GM/r 3 ) 1 is the Keplerian angular velocity. Typically the contribution of 
self-gravity gives Q — Qk = O^fl^r^ 1 ). As a check, we evaluated the contribution of the pressure 
term for the models we consider here and verified that it was relatively small. 

We could also add the time average potential due to a companion with small mass, which orbits 
in the disk plane, into ^* ex t and incorporate it in the equilibrium. This would give a contribution 
of order u z , being the orbital precession frequency of a free disk particle, to $7 — Qk- However, 
we shall limit ourselves to companion masses small enough so that \u z \ <C (HQr^ 1 ) and its effect 
may in general be considered to be small in comparison to that due to self-gravity. But, in a 
formal sense, we shall consider that typically \u z \ 3> (H 2 £lr~ 2 ), so that it may be retained while 
frequencies of order H 2 £2r~ 2 may be neglected. 
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The warping and precession of a disk at small inclination may be considered as a small amplitude 
perturbation about such a state incorporating a companion (e.g. Toomre 1983; Sparke 1984) and 
see below. 

In principle, the surface density distribution, S = ff^ pdz, is determined by the viscous evolution 
of the disk. In the absence of a deterministic prescription for the effective viscosity, for illustrative 
purposes we arbitrarily choose S(r) = So (R/r — 1) , So being a constant. Alternative functional 
forms describing the way S tapers to zero have also been considered below. 

2.4. Linear perturbations 

In the limit where either any perturbing/forcing potential is small compared to that due to the 
central mass and the amplitude of any free bending modes is small, the disk response can be 
described in a linear analysis. We here follow the standard procedure of vertical averaging (see, 
for example, Hunter & Toomre 1969, Sparke 1984, Papaloizou & Lin 1995). This should be valid 
when the radial wavelength of the disturbance is significantly longer than the disk thickness as is 
the case for the perturbations considered here. For linear warps, the dependence on <p and t may 
be taken to be through a factor ex.p[i(m(p + at)]. Henceforth this factor will be taken as read. For 
the time being we shall leave it as a multiplicative factor for all perturbations but will remove it 
at a later stage. The mode frequency, a, can be taken to be zero for a secular response viewed 
in an appropriate rotating reference frame and we shall consider this case below. Alternatively, 
the disturbance may appear with a small non zero frequency or pattern speed when viewed in an 
inertial frame. In the latter case, the disk will appear to precess if the perturbation has a global 
form. We limit ourselves to considering frequencies significantly less than the inverse density wave 
crossing time of the disk (characteristically 0(ffOr -1 )) which is a necessary condition if the disk 
is to exhibit global precession (Hofner & Sparke 1994, Papaloizou & Terquem 1995). Then we 
show below that, for modest warping of the disk, these two descriptions are equivalent. Finally, 
we remark that the azimuthal mode number m = 1 throughout. 

We denote the Lagrangian displacement as £ = (£ r , £ z ). Here we shall not assume at the outset 
that the horizontal components of the displacement are zero, as it was shown in Papaloizou & 
Lin (1995) that these can produce significant effects in a self-gravitating but near Keplerian disk 
of the type considered here. 

The inclination (vertical displacement /radius) associated with the tilt of the disk is £ z /r. This is 
in general assumed to be a slowly varying function of both r and z such that its variation with z 
may be neglected. The applicability of this vertical averaging approximation is supported by the 
numerical calculations of Papaloizou & Lin (1995) which allowed for the possibility of variation of 
the vertical displacement associated with modes of the type we consider with z. 

We first assume that we are working in an inertial frame in which the unperturbed disk appears 
steady. However, it may subsequently be convenient to adopt a uniformly rotating frame so as to 
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remove the rigid body precession associated with the motion of the disk, in which case the effects 
due to the coriolis force need to be added. 

We begin by writing down the perturbed form of the vertical component of the equation of motion 
applicable to a gaseous disk with a barotropic equation of state (see Hunter & Toomre 1969) 



D% = _d^ _ d{P'/p) _ 
Dt 2 dz dz dz 



(10) 



where perturbations to quantities are denoted with a prime and the convective time derivative 
D /Dt is here equivalent to multiplication by i{a + fi). 

The potential perturbation, *&' ext , is taken to be the secular contribution due to the companion in 
an inclined circular orbit (see § |2.2| ) where this has not been included in the equilibrium. Then 
we calculate the disk response taking a = 0, but will need to transform to a rotating frame (see 
below). If the companion orbits in the initial symmetry plane of the disk and its effect included in 
the equilibrium, then we look for free bending modes with small frequency 0(u z ), the amplitude 
of which will give the inclination of the precessing disk. These two approaches are equivalent for 
small inclinations and modest warping. 



For a polytropic equation of state we also have P' = p'c 2 s . Using the above, equation (10) may be 
written 



r^2, d*'r d (p' c2 s/p) dKxt 
dz dz dz 

The perturbation to the gravitational potential of the disk is given by the Poisson integral 



V |r — r' 



*» = -a (12) 



while the perturbed continuity equation gives 



p' = -v ■ G*) . (13) 

We then find after an integration by parts, assuming that the disk density vanishes at its 
boundaries, that 

m' G = —G f ?(rWHr-r-HV. (14) 
Jv Ir-r'l 



In what follows, we find it convenient to consider the contributions to p' and ^f' G arising from the 
horizontal and vertical components of £ separately. The contributions arising from the vertical 
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component we denote with a subscript v and call them vertical contributions. The contributions 
arising from the horizontal components we denote with a subscript h and call them horizontal 
contributions. Thus 



p' v 



dz 



^?'q v and ^' Gh denote the potential arising from p' v , and p' h respectively. Thus 



(15) 
(16) 



Gv 



Gh 



-G 



-G 



P (rQ6(rQ (z-z')dV 
\ |r-r'| 3 
p (r') £ h (r') • (r - r') dV 
Ir-r'l 3 



v 



(17) 
(18) 



The vertical and horizontal contributions give rise to additive contributions, and F' h respectively, 
to the vertical component of the equation of motion (11). We consider these in turn below. 



2.4-1- The tilt equation 

Having evaluated the vertical and horizontal contributions we may construct a tilt equation 
governing the inclination g derived from vertical integration of equation (11) in the form 

/co roo F)\Ul 

p{F' h + F' v )dz+ 9-^dz. (19) 
-oo J —co OZ 

2-4-2. The vertical contribution 
From equation (|Tl|), this is seen to give rise to the vertical deceleration 

oz Oz 



To deal with this, we differentiate the vertical contribution to the Poisson integral ( |17|) with 
respect to z, differentiate the Poisson integral (|3|) twice with respect to z and combine with £ z to 
obtain in the limit of an arbitrary thin disk in its midplane 

— + iz ^^- G J A b~T^ < W 
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where the two dimensional integral is over the surface of the disk. 
Proceeding further, we remark here that we shall eventually require 



pF' v dz 



S(r') [L(r')-L (r)]dV 



dz- 



oo Q 

P 



] d{ptj z 
dz \ p dz 



j dz. 

(22) 



Here we have made use of ( |2l| ) and (|15|). Equation (22) may be further simplified by using 
the assumption that £ 2 is independent of z, together with use of the fact that the equilibrium 
quantities satisfy hydrostatic equilibrium (@) in the vertical direction, with the result that 



pKdz _sf.ggjj._fl/SMSC) fc M - f. (r)l _V 
0£ |r — r'l 



Here and below, we make the assumption that the second derivatives of the external potential 
are slowly varying through the disk thickness so they can be evaluated in the midplane. This 
completes the reduction of the vertical contribution to the force per unit area which now depends 
only on the vertical displacement. 



2.4-3. The horizontal contribution 
From equation (|Tl|), this is seen to give rise to the horizontal deceleration 



F , = s^ + sjMiA (24) 

dz oz 

Noting again that we need to perform a vertical average, we multiply by p and integrate over 
the disk thickness to obtain after an integration by parts of the second term and use of vertical 
hydrostatic equilibrium 



L pFhdz = L v— +ph ^z~) dz + loo Ph — dz - (25) 

We pause for a moment to consider the magnitudes of the terms in equation (|2q). Referring 
back firstly to equation (|23|), we comment that the first term on the right hand side there is of 
order f_^__£ z while the second, being produced by the self-gravity of the disk, is characteristically 
smaller by a factor H/r for Toomre parameter Q of order unity and global perturbations as 
considered here. Allowing for the possibility that the horizontal displacement could become as 
large in magnitude as the vertical one, the first two terms in equation (p_|), being produced by 
the disk self-gravity, are each also of potential magnitude __^Si?^ 2 /r. However, we show in the 
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Appendix that, for global perturbations, the combination is in fact smaller by an additional factor 
of order H/r, and so can be neglected in our ordering scheme. Thus only the final term remains 
in equation (pl|). This is best dealt with by considering the perturbed equations of motion in 
Eulerian form. These can be written (remembering m = 1) 



i{(7 + n)v! r -2Qv' ip = 

K 2 

i(a + n)v' ip + —v' r = 

i(a + 9)v' z = 



dW 
dr ' 

iW 
r 

dW 

dz 



(26) 
(27) 
(28) 



Here the velocity perturbation is v' = (v' r1 v'v' z ), W = P'/p+ and k 2 = 20r l d(r 2 £l) /dr is 
the square of the epicyclic frequency. 

Given that we seek an inclination £ z /r = v' z / [ir (a + f2)] that is independent of z, we may write 
W = (a + n) 2 z£ z . Equations © and © may then be used to solve for the horizontal velocity 
perturbations in terms of £ 2 . These can then be used to find p' h from the perturbed continuity 
equation (|l3|). An analagous procedure has been carried out in Papaloizou & Lin (1995) and 
Papaloizou & Terquem (1995). As we are here interested in the low frequency limit, it is consistent 
with our ordering scheme to perform this procedure adopting a = 0. We thus obtain 



, = -xiz d(tj z r 2 n 2 ) 

r 2 (K 2 -n 2 ) Or ' K ' 

k 2 z d(g z r 2 n 2 ) _ z^n 

V * 2r 2 Q(K 2 -n 2 ) Or r ' 1 j 



Finally we obtain from (|13| ) 

Ph in {P h) OrV2 dr \r*/ 2 (k 2 - tl 2 ) dr )' 1 j 

In forming the above, consistently with the neglect of frequencies of order H 2 £l,K r ~ 2 , we have 
neglected the contribution of the last term in (EH) to (|3l|) . 



Using equation (p5|), we find the horizontal contribution to the vertically integrated vertical 
component of the equation of motion 

w, _ d 2 ^ ext 1 9/ fjZl d(£ z r 2 n 2 ) \ 
3 9 h ~ dz 2 UrV 2 [ r3/2 («2 _ n2 ) Q- r )' ^ 
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where we have approximated d^ ex t/dz = (d 2 ^ ex t/dz 2 )z, the second derivative being evaluated on 
the midplane, and \i = pz 2 dz. 



2.4-4- Reduction of the tilt equation 



Having evaluated the vertical and horizontal contributions, we construct the tilt equation (|I 
which now appears as an equation for £ z alone. 

Using equations (p3[) and (|32[), we obtain 



E(r)S(r') [L(r')-L (r)] d 2 r' 



n2 d ( mO| afe/r) 

<9r \ k 2 — f2 2 dr 



e.ct 



6*2 



(33) 



Here, consistently with the neglect of frequencies of order H 2 £li< r ~ 2 > we have made the 
replacements d 2 ^ ex t/dz 2 = $l 2 = Q 2 K in ( |32l) everywhere apart from in the denominator k 2 — $7 2 . 
Also in the last term, the derivative is evaluated on the midplane. 

Rather than retaining £ z as the quantity to be determined, we adopt g = i (£ 2 /r) exp (— iip) . 
At zero frequency, its modulus is the ratio of vertical to azimuthal velocity and also the disk 
inclination. Remembering that the azimuthal dependence of £ z is through a factor exp (i<p) , 
the replacement with g which is accordingly independent of <p will remove this factor from 
equation (|33"|). Carrying out this reduction we obtain 



G 



d 2 ^ 



ext 



r dr 



8z 2 



n 2 K ~ L{9) + r n 2 K dz ■ 



where the operator C is defined through 



(34) 



g r R 



C (9) = -1^2- I K (n r ') S(r)S(r') [g(r) - g[r')) r' 2 r 2 dr' + 



dr \ k 2 — Q, 2 dr 



(35) 



Here the kernel K is given by 



K(r, r') 



277 



cos <M<!> 



[r 2 + r' 2 — 2rr' cos <3?] 



3/2 ' 



(36) 



and R; is the innermost disk radius. 
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We also note that after using (]?]), for the low frequencies we consider, to the accuracy we are 
working 

n2 id^> G d 2 v ext id^ ext d 2 v ext on 

^ 5 X~2~ = 5 Wl~ = Ill K^z, 

r or oz z r or oz z 

where again u z (r) is the orbital precession frequency due to the companion if it orbits in the 
symmetry plane and is included in the equilibrium external potential so that ^' ext = 0. If instead 
the effect of a companion in an inclined orbit is included only through ^' ext , we have oj z = 0. 

Neglecting a 2 , the inviscid tilt equation is then written in the convenient form: 



n K {9) + r n 2 dz 



Aff)+ „2 (37) 



We comment that the above equation describes the global response of a thin near Keplerian 
self-gravitating disk with Toomre parameter Q of order unity. The largest frequencies neglected 
are of order H 2 Q,k r ~ 2 - It can also be used to describe low frequency bending modes when the 
forcing term is set to zero. There are two identifiable contributions to the operator C defined 
by (H). The first comes from self-gravity and leads to the description of warps given by Hunter & 
Toomre (1969), Sparke (1984), Sparke & Casertano (1988), Kuijken (1991) and others. The 
second term in C can be identified as arising from pressure and it can lead to comparable effects 
to those due to self-gravity in a near Keplerian disk, in contrast to a galactic disk, because of the 
near Lindblad resonance that arises in the former case because k — Q = 0(H£lr~ l ). The form of 
this term is the zero frequency global counterpart of that occurring in the low frequency WKB 
dispersion relations given by Papaloizou &: Lin (1995) and Masset &: Tagger (1996). 

Another useful property noted by Sparke & Casertano (1988) in the pure self-gravity problem is 
that, for a finite isolated disk, C is self-adjoint, that is for any pair (51,52) that satisfy appropriate 
regularity conditions, including zero derivative at a Lindblad resonance if the pressure term is 
included: 



/ g* 2 C( gi )dr= / glC(g 2 )dr) . (38) 



Furthermore, when uj z = and there is no forcing, equation ( p7| ) has the solution g = constant 
which coresponds to a rigid tilt in a pure Keplerian potential. Modification of this rigid tilt mode 
by a companion results in a low frequency modified tilt mode that describes global warping and 
precession of the disk. 
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2.4-5. Disk precession 

i) Low frequency bending modes 

We first consider the case of a disk precessing at low inclination about the angular momentum 
vector associated with the companion orbit. This can be viewed as the excitation of a small 
amplitude bending mode with the time averaged potential of the companion being included in the 
equilibrium state. Implicit in this treatment is the notion that the companion orbit is fixed so 
its angular momentum content must be presumed to be significantly larger than that of the disk. 
Such modes are governed by the eigenvalue problem 



2(. + ^)S g =£to) (39) 

When uj z = 0, there is a solution corresponding to the rigid tilt mode with g = go = constant = 8, 
where for convenience we take 5 to be real and therefore equal to the inclination and a = 0. For 
small uj z , there is a solution with g close to go an d <r non zero but small in magnitude. From first 
order perturbation theory, remembering C(go) = 0, we obtain the forced response problem for g: 



_ (40) 

\l K 

Multiplying by go, integrating over the disk and using the self-adjoint property of C, we find the 
integr ability condition which determines a: 



a f^ dr = -f R ^ dr . (41) 
JRi JRi Mr 

Remembering that go is constant and the form of Kepler's law, this may be written 

<t = - f 2-KrY,{r)oj z {r)j{r)5dr I f 2irrT l (r)j(r)5dr , (42) 
J Ri / J Ri 



where j(r) = r 2 VLx denotes the specific angular momentum of the orbiting disk material. 

This relation has a simple interpretation, as the gyroscope equation, that the magnitude of the 
disk precession frequency —a should equal the component of the applied torque perpendicular 
to the companion angular momentum axis and the disk angular momentum axis divided by the 
magnitude of the component of the disk angular momentum vector perpendicular to the companion 
angular momentum axis (Kuijken 1991). This relation is applicable to secular precession of free 
particle orbits and thus determines oo z , so u> z j(r)8 = —T(r), where T(r) is the external torque 
component per unit mass acting on a disk circular orbit at radius r. 
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Once the precession frequency has been found, the forced response problem can be solved to find 
the degree of warping. For the idea of global precession to be valid this should be small. This 
condition is roughly equivalent to the requirement that the disk wave crossing time be short 
compared to the precession period. 



Another usefull relation that can be obtained from fl39| ) is found by multiplying by g* and 
integrating over the disk, giving 

2a [ ^-^-dr + 2 f ^ g 9 dr = [ g*C(g)dr. (43) 



The eigenfrequency a, if real, corresponds to a simple precession. However, if it were complex, 
the imaginary part would give the rate of decay of the inclination towards the symmetry plane. 
From ( f43| ) we obtain for this decay rate 

Im{a) = Im {J g*C(g)dr^j j (2 J -fi-^drj . (44) 

Because of the self-adjoint property of the problem considered up to now, the decay rate is zero. 
However, if contributions from viscous or other non conservative forces are added, a non zero 
decay rate may result (see below). 

ii ) External forcing at finite inclination 

In this case, we consider the zero frequency response of a disk forced by the secular potential due 
to a companion in an inclined orbit given by (^) so that uj z = a = 0, and the governing equation is 

c{9) = -rt5"af ■ (45) 

But because the unforced problem has g = constant as a solution, the forced problem will not 
in general have a solution. Physically, this is because there is an unbalanced torque due to the 
companion which produces disk precession. To deal with this, we suppose the disk appears steady 
in a frame precessing about the total angular momentum axis with frequency u„. At first, let us 
suppose this coincides with the companion orbital angular momentum vector as in the preceeding 
section. Then we must add the vertical component of the coriolis force as a perturbing force along 
with that due to the external potential. This amounts to an additional perturbing force per unit 
mass being the real part of — 2iuj p sin(<5)rOjf exp (i<p) . Incorporating this into equation ( |40l) (i.e. 
adding it to —d^>' ext /dz) results in this now becoming 

C{9) = d?£ ( 2 ^Ksm6- i^4) . (46) 
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The above problem is very similar to, and in fact equivalent to, ( |40| ) for small S. Just as in that 
case, we have a forced response problem to solve for g and we find the integrability condition that 




determines u) p to be 

— dr. (47) 

Ri Mk 

This equation can also be interpreted as a gyroscope equation giving the precesion frequency as 
the ratio of appropriate disk torque and angular momentum components. Note here that the 
perturbing external potential is complex (see |6|) ensuring that the calculated precession frequency 
is real and that the factor of two in (f47|) accounts for an azimuthal average. 



After finding uj p from (47), equation (46) may be integrated to give g, to within the addition 
of an arbitrary constant inclination, which may be eliminated by choosing the coordinate system 
so that g = at the disk inner boundary. Then if g is small elsewhere, the disk approximately 
precesses like a rigid body. 



We have mentioned in § [L3] that the disk is not expected to have settled to equilibrium with zero 
inclination to the orbital plane. The choice of zero inclination perturbation at the inner boundary 
is then justified. 

From (||), we see that uj p oc cos 5 when D 3> R such that the tidal potential becomes quadrupolar. 
The same dependance is found for a disk subject to the potential of a misaligned flattened halo, 
5 being in that case the angle between the disk plane and the symmetry plane of the potential 
(Kuijken 1991). 

Up to now we have assumed that the disk angular momentum content can be neglected in 
comparison to that of the orbiting companion. Then the disk precesses about the companion 
angular momentum axis. For finite disk angular momentum content, the disk and companion 
both precess about the conserved total angular momentum vector. This can be taken account of 
by noting that the precession frequency uj p as calculated above should be replaced by ui p J or b/J, 
where J or b is the companion's orbital angular momentum and J is the total angular momentum. 



2.4-6. The effect of viscosity 

In principle viscosity can act on the disk so as to change its inclination with respect to the orbital 
plane (see § |1.3| ), and also by damping the horizontal velocities that can be amplified by a close 
Lindblad resonance (see equations ^ and [3C]) . 



We have already observed in § 1.3 that, as long as it is globally warped, the timescale on which 
the disk inclination is expected to change is the global viscous timescale, namely Q~ 1 (R/ H) 2 /a, 
where a is the standard Shakura & Sunyaev (1973) parameter. Since a < 1, evolution on this 



-21 - 



timescale is long and has been neglected in comparison to the other processes considered in the 
above analysis, where frequencies on the order of H 2 £lr~ 2 have been neglected. Thus, in the 
ordering scheme we have adopted, changes in the disk inclination with respect to the orbital plane 
should be negligible. 

In a near Keplerian disk, viscosity also provides a damping effect on the near resonantly produced 
horizontal velocity perturbations. We follow the discussion given by Papaloizou & Lin (1994, 
1995). The near resonant denominator acting in equations (p9|) and (|30| ) is 1 — k 2 /£1 2 ~ 0(H/r) 
when Q ~ 1. The effect of an appropriately defined a viscosity, at the order we are working, is to 
replace this denominator by (1 — ia) 2 — k 2 /Q 2 . Demianski & Ivanov (1997) have recently given the 
relativistic generalization. Thus viscosity can only be neglected when a <C H/r, which is not the 
case in some models of NGC 4258 (see below). 

As indicated above, we take viscosity into account by changing the denominator k 2 — O 2 of the 



last term of the operator £ (equation 5^) to k 2 — Q 2 (l 



ia) 



We note that when a 3> H/r, 



the pressure term in C is small compared to that due to self-gravity, so that the warp is mainly 
controlled by self-gravity. This is illustrated below in the case of models for NGC 4258. 

We comment that incorporation of viscosity as described above does not affect the discussion 



leading to the determination of the precession frequency through equation (47). This is consistent 
with the idea that as long as the disk is only mildly warped, the internal viscous forces do not 
enable a net torque to be exerted on the disk, and so do not modify the precession rate or change 
the disk inclination except possibly on a long timescale. 



In the case of small inclination we may use (44) to find the rate of decay of the inclination after 
incorporating the effect of viscosity as outlined above. For this we then obtain 



Im(a) 



[k 2 - Q 2 (1 - a 2 )] z + 4a 2 ft 4 



dg_ 
dr 



dr 



Zg*g 

n K 



dr 



(48) 



At a given location in the disk, the decay rate is a maximum for a ~ H/r. For a > H/r, the 
decay rate is found to be of order Q,(H/r) 2 a [5g/(ag)] 2 , where 5g represents the total change in 
inclination or the total degree of warping in the system. Thus the decay rate is small if the total 
warp is small. 



3. Application to NGC 4258 

A warped disk model was originally proposed for NGC 4258 to account for the origin of the 
megamasers (Neufield & Maloney 1995). In this scenario, the warp provides a favorable condition 
for parts of the disk to be illuminated by the X-ray source within the nucleus. In the regions 
between 0.13 — 0.26 pc, the surface layers of the disk are heated to temperatures in the range 
300 K < T < 8,000 K at which water maser production may occur (Collison & Watson 1995). 
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Neufield & Maloney (1995) supposed that conditions are less favourable for reprocessing in those 
regions of the disk interior to 0.13 pc and exterior to 0.26 pc. 

With respect to the systemic velocity of NGC 4258, the flux of the red-shifted masers is an 
order of magnitude larger than that of the blue-shifted masers. Herrnstein et al. (1996) suggested 
that the heated gas above the concave side of a warped disk leads to greater extinction than that 
above the convex side. With the appropriate orientation, thermal absorption is dominant for the 
blue-shifted masers. 



3.1. Disk Parameters 

For the total luminosity of L = 4 x 10 41 erg/s attributed to this object, an efficiency e (here 
taken to be 0.1) of rest mass conversion implies a mass accretion rate through the system of 

M=4=6x 10~ 5 M Q /yr. (49) 



ec 



Note that the need for reprocessing is evident as according to a steady state optically thick 
disk model (see Pringle 1981) the effective temperature would be given by 

T eff = ll(r 18 r - 75 K, 

where rig is the radius r in units of 10 18 cm. For comparison, the effective temperature, assuming 
that the central luminosity is radiated over a radial scale rig, would be 

T eff = 154(n 8 )- - 5 K. 

Thus, in the most region of the disk, the irradiation is the dominant source of heating and T e ff is 
in the range appropriate for water maser emission. 

The disk semi-thickness H is given by 

^ ~ ^ ~ 1.5 x 10- 3 p^V /2 • (50) 

Here M7.6 = M/(4 x 10 7 M ), M is the central black hole's mass, c s is the sound speed, 
T200 = T c /200 K, T c is the midplane temperature, and Q is the orbital rotational frequency in the 
disk. From the observed upper limit H/R < 2.5 x 10 -3 , we find T c ~ T e ff due to irradiation and 
T 2 oori8/M 7 .6 < 2.8. 

Applying the conventional a prescription for an effective turbulent viscosity (Shakura & 
Sunyaev 1973), the mass transfer rate in a steady-state disk is 
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M = 3iraH 2 nZ = 0.59a— ( ^ ur "''° I M /yr, (51) 

M \ r 18 J 

where we have replaced 7r£r 2 , with £ being the surface density, by Mp, the disk mass interior to 
radius r. For NGC 4258, we find from (|l|) and © that 

4 M /T 9 2 nn M 76 \~ 1/2 . , 

a - 10" 4 200 7,b . 52) 

M D \ ns / 

The importance of self-gravity is measured by the Toomre Q parameter, Q ~ (HM)/(rMo), such 
that Q > 1 is required for stability to axisymmetric ring modes. From (|50| ) we find that the disk 
becomes gravitationally unstable in this way if M D /M > 1.5 x 1(T 3 . From (||) we find a ~ 0.1 
for Mjy/M = 10~ 3 . A smaller value of a would imply larger S, small Q value and a more unstable 
disk. 

Based on the above estimates, we consider a model which is marginally self-gravitating with Q ~ 1 
and a ~ 0.1 such that M/) ~ 3 x 10 M Q . We note that Maoz (1995) has suggested that the radial 
intervals between successive maser sources may be related to the wavelength of transient spiral 
structure in a marginally self-gravitating disk. The value of a we have assumed is comparable to 
that expected to occur in a marginally self gravitating disk as a result of non-axisymmetric waves 
(Laughlin & Rozyczka 1996). 

The inner part of the disk of NGC 4258 has been modelled by Lasota et al. (1996) as being 
advection-dominated. In this case there would be greater mass flow rate through the disk than 
indicated by the observed luminosity. In principle this situation could be obtained by taking a 
somewhat thicker more massive disk together with a larger value of a. As indicated below, our 
results may be scaled to apply to such a case. 

Here, we consider the response of a self-gravitating disk, with parameters like those discussed above, 
to an orbiting binary companion in an inclined circular orbit. This interaction produces warping 
and precession of the disk. We note that angular momentum transfer between the companion and 
disk may produce truncation and gap formation in the disk (Lin & Papaloizou 1993), terminating 
the radial progression of maser sources. 

Pringle (1996) has indicated that forces due to radiation pressure may be important in producing 
warps in disks. Such forces could be included in the formalism presented here. However, we shall 
limit ourselves in this application to discussing the effects arising from a companion only. 

In our model, the companion then needs to have a mass comparable to that of the disk and an 
orbital radius somewhat larger than the size of the disk. Black hole binary systems in which the 
companions can have up to comparable mass and separations on a scale of 10 18 cm have been 
considered by Begelman, Blandford & Rees (1980), and postulated to account for jet behavior 
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in extragalactic sources by Kaastra & Roos (1992) and Roos, Kaastra & Hummel (1993). The 
dense cores of globular clusters may be an alternative candidate for such a companion. From 
the negative intrinsic period derivative of millisecond pulsars, Phinney (1992) derived a mass 
> 4.5 x 10 3 Mq which is concentrated within 0.05 pc from the core of M15. It is possible that 
such a core would not only be adequate to provide the necessary perturbation for the observed 
warp but also survive the tidal disruption by the central black hole. 

3.2. Numerical Results 

For illustrative purposes, we present a numerical model to show that the observed properties of the 
warped disk in NGC 4258 can be regulated by the perturbation on a marginally self-gravitating 
disk due to a low-mass companion on an inclined orbit. Since we include the secular effect of the 
companion in an orbit with finite inclination through the perturbing potential ^' ext , we consider 
equation ( |34"1) with u z = (i.e. the left hand side is zero). 

In our numerical calculations, we solve equation (|46| ) by dividing the interval [Ri,R] into a grid 
( r j)j=i n usrn § a constant spacing Ar in the radial direction and such that n = R% + Ar/2 and 
r n = R — Ar/2. We approximate equation ( |34] ) at r = rj, j = 1 . . . n as a system of n linear 
equations for the n quantities gj = g(fj),j = 1 . . . n (typically n = 200-2000 in our calculations). 
The fact that gj = constant is a solution gives a solubility condition analagous to ( f47|) from which 
the precession frequency may be found. The perturbing potential can then be modified so as 
to include the effects of the Coriolis force. The specification g\ = 0, corresponding to zero disk 
inclination at the inner boundary, then enables the solution to be completed. 

For computational convenience, we normalize the units such that M = 1, R = 1 and 
Qk(R) = 1- The observations of Miyoshi et al. (1995) indicate that Md(R) of the disk is less 
than 10 _2 M. Marginally self-gravitating disks are then such that the disk semi-thickness is less 
than one percent of the radius. For illustrative purposes we have chosen a model for which 
S = £o(i?/r — 1), with the constant So being chosen such that Md(R) = 1.5 x 10 _3 M, a = 0.1, 
and the maximum value of HJt^ (-f^/T)maxi 

to be ~ 10 3 (see discussion in § [j]). 

The equilibrium structure of the disk is calculated adopting a polytropic index 1.5. We plot 
the Toomre parameter Q = kc s /(7tGX), with c s being evaluated in the midplane for our model 
disk with R-i = O.li? in Figure |]. The minimum value of Q is approximately unity making the disk 
marginally self-gravitating. 

For the companion providing the perturbation, we adopt an orbital inclination with respect 
to the disk of 5 = ir/4. Results for other values of 5 can be obtained by noting that to a reasonable 
approximation, the precession frequency is proportional to cos 5 and g is proportional to sin 25. 
We also take the companion's mass to be such that M p = Md(R) and the orbital radius to be 
D = 1.5R. 
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From equation (p3|)) we note that g is essentially independent of a scaling which reduces 
Md,M p , and H by the same factor. This scaling relation holds because Q is unaltered. However, 
for large a > H/r, the small amount of twist induced by the viscosity is reduced by the same 
factor. Other things being equal, g scales with M p . Results were found to be insensitive to the 
form of the surface density taper at the outer boundary. 

The induced precession frequency in our model is uj p = —7.9 x 10 -5 (J/ J or b)Q(R). We are 
able to reproduce a warp with a relatively large amplitude using a perturbing potential which is 
much smaller than that due to the central mass, justifying the assumption of a small perturbation. 
An approximate estimate for \g\ may be obtained by estimating the magnitudes of terms in (j6|) 
and (|34|) as 

\q ~ 7—; — 5- sm2o . (53) 

This amplitude can be interpreted as being the product of the companion induced precession 
frequency and the time required for a bending wave to propagate through the disk. When this 
product is small, there is good communication between the different parts of the disk, and the 
warp is small. In the thin, self-gravitating disks we have considered here, wave propagation is 
mainly regulated by self-gravity which acts to preserve the disk's intrinsic spin vector. 

In Figure ||.a, we plot both the real part of g and the imaginary part of —g as functions of 
r. These quantities represent the inclination of the disk as seen edge on, viewed looking down 
the y and x-axis, where <p = tt/2 and respectively. We have also calculated g with just the 
self-gravitating term in (|35|). As expected (see § 2.4.6| ), the real part of g so obtained is hardly 



distinguishable from what we get when pressure and viscosity with a = 0.1 are included. For 
comparison, Figure |2[b shows the real part of g obtained including both self-gravity and pressure 
terms but with very small a = 10~ 5 . When a = the solution looks almost exactly the same. 
The case a = 0.1 has been displayed again on the same plot for comparison. We observe that 
when a is small, the response oscillates around that associated with the larger a. This is due to 
the fact that when both self-gravity and pressure are present, both a slow and a fast bending 
wave can propagate in the disk (Papaloizou &: Lin 1995). When a is large, the pressure term 
in (|35| ) is small, so that only the fast wave corresponding to the purely self-gravitating case can 
propagate. In Figure §.b we see the superposition of the fast and slow waves, the former having a 
larger wavelength than the latter. We have checked that the slow wave disappears once a becomes 
larger than H/r. We note that the precession frequency is the same in all cases. Hereafter, unless 
otherwise specified, we shall refer to the case where self-gravity, pressure and viscosity with 
a = 0.1 are included. 

In Figure ||, we plot a three dimensional view of the warped disk. For clarity, the vertical 
scale has been magnified. 

The observations indicate the presence of a warp corresponding to \g\ reaching values of about 
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a few tenths at the disk edge. The results presented here reproduce a warp of the required 
magnitude. Since the outer region of the disk is mainly regulated by the external perturbation, 
the functional dependence of g on r in the outer part of the disk is insensitive to the choice of Ri. 
But the total range of \g\ increases as Ri decreases. This dependence is due to the fact that in our 
model the wave crossing time across the disk decreases as Ri increases. In the limit of small Ri, it 
becomes difficult for the outer parts of the disk to "communicate" with the inner parts such that 
the two regions become essentially disconnected. 

We now consider the twist (the azimuth of the maximum of the real part of g at a given radius 
as a function of radius). This azimuth is given by 2tt — rj(r), with r)(r) being the argument of g. 
Figure |I| shows r/(r) as a function of r. 

The total variation in r\ is about a few degrees and most of this variation occurs in the inner 
regions. A positive gradient dr)(r)/dr implies a negative gradient in the azimuth and a regression 
of the line of nodes which corresponds to a trailing twist. This trailing twist is produced by 
the second term of ( pop when viscosity is added (see § 2.4. 6| ). In the outer part of the disk, the 
magnitude of this term is a factor H/r smaller than the other terms. If it is completely neglected, 
as noted above, \g\ is hardly changed, indicating that its only significant contribution is to produce 
the small twist. The twist becomes larger in the inner part of the disk because the perturbation 
amplitude of the external forcing is weakened relative to the effect of inertial response and viscous 
stress which are contained in the second term of (|35|). The value of rj at small r would then be 
smaller if Q were not as large there as it is in our model. 

Accurate determination of the twist can provide useful data on the magnitude of a. According 
to the results in Figure Q the magnitude of rdr]/dr ~ 1 in the outer part of our model disk 
(with a = 0.1) and an order of magnitude larger in the inner part. In the analysis of NGC 4258, 
Herrnstein et al. (1996) constructed models to fit the observed position and velocity of the maser 
sources. Some of their models include variations, along the radial direction, in both the inclination 
and the projected position angles of the disk's spin axis. For example, their model 3 is consistent 

o 

with a small trailing twist in which \rdrj/dr\ ~ 2 , although models without any twist fit the data 
equally well. We note that the warp amplitude is on the order of a tenth of the radius throughout 
this particular model. 



4. Discussion 

We have derived the tilt equation governing the inclination of a near Keplerian disk subject 
to global warping (m = 1) perturbations taking into account self-gravity and pressure. We have 
considered the case where the perturbation is due to a companion on an inclined orbit, but the 
analysis we have presented in this paper is valid for any m = 1 perturbing potential proportional 
to the coordinate z along the disk rotation axis. All effects that can cause evolution on the fastest 
possible timescale on the order of HQr~ l have been included. We also considered the effects of 
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viscosity, showing that pressure effects diminish in importance relative to those due to self-gravity 
once the viscosity parameter a > H/r. The disk response and precession frequency have been 
derived using two different approaches that we have shown to be equivalent in the small inclination 
limit. We considered i) the excitation of small amplitude low-frequency bending modes about the 
orbital plane of the companion. Here, the time averaged potential of the companion is included 
in the equilibrium state, ii) the zero frequency response of a disk forced by the secular potential 
due to a companion in an orbit with finite inclination. It was found that when the precession 
frequency is small compared to the rate of wave propagation across the disk, the disk can undergo 
quasi rigid body precession with a small warp. Then the timescale for the inclination to decrease 
is found to be long compared to other timescales in the problem. 

For an illustrative application, we have explored the possibility that the recently observed 
warped disk in the active galaxy NGC 4258 is produced by a binary companion. Our results 
indicate that it can be produced by a companion with a comparable mass to that of the observed 
disk. For a disk mass of 10 -3 M = 4 x 10 4 M , the required mass for the companion is comparable 
to that contained in the cores of the densest globular clusters. Either these cores or black holes 
with comparable masses can survive the tidal disruption of the central black hole. The warp 
remains of modest magnitude as long as the companion induced precession time is everywhere 
significantly longer than the time required for a density wave to propagate through the disk. 

If the disk and companion's mass are comparable, their tidal interaction would lead to the 
tidal truncation of the outer disk edge and the companion's orbital evolution on the viscous 
timescale of the disk. The exchange of angular momentum between the companion and the disk 
occurs through the excitation of density waves in the disk. For the parameters we adopted, the 
viscous evolution timescale of the disk is ~ 8 x 10 s yr giving a significant lifetime in the current 
state. 

Subject to the effect of dynamical friction, the companion may also undergo orbital decay as 
it interacts with the field stars in the nucleus of NGC 4258. For illustrative purposes, we adopt 
a model in which the density of the galactic core varies like r~ n , r being the distance from the 
galactic centre. Assuming that the velocity of the companion stays Keplerian as it spirals down to 
the centre, the timescale td on which the angular momentum L of the companion changes is then 
(Binney & Tremaine 1987) 

L C 1 M M 
td ~ dL/dt ~ 4^uM pPo D 3 

where po is the stellar density at the distance r = D from the centre. C is a number which depends 
on various parameters of the galactic core (see Binney & Tremaine 1987 for more details). It is 
typically on the order of unity, and we then set C = 1. In our model, D ~ 0.3 pc. If we take for po 
the value estimated for the galaxy M32, which has been resolved down to subparsec scale, we have 
Po between 10 6 (Lauer et al. 1992) and 10 7 M & /pc 3 . This gives td between 10 7 and 10 8 yr, and 
thus the companion loses angular momentum because of dynamical friction faster than it gains it 
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by interaction with the disk. For td to be at least on the order of the viscous timescale of the disk, 
we need po < 1.5 x 10 5 M & /pc 3 . This is not unreasonable, since we might expect the stars in the 
centre of NGC 4258 to be cleared out by the presence of the companion. 

The fact that the companion spirals in with a lifetime comparable or smaller than the viscous 
timescale of the disk provides an additional argument why we do not expect the disk to have 
settled to equilibrium with zero inclination to the orbital plane (see discussion in § |l]). 

The authors wish to thank L.J. Greenhill for useful conversations. This work is supported by 
PPARC through grant GR/H/09454 and by NASA through grant NAG 53059. C.T. acknowledges 
support for this work by the Center for Star Formation Studies at NASA/ Ames Research Center 
and the University of California at Berkeley and Santa-Cruz, and by the EU. 

A. APPENDIX 

We here show that the magnitude of the combination of the two terms in equation ( p5|) given by 

for global perturbations under conditions for which the Toomre parameter Q ~ 1, is of order 
Q 2 K Y l H 2 £ >z /r 2 for characteristic values of H and r and so can be neglected in our ordering scheme. 
Of course, because these terms derive from the disk self-gravity, for Q > 1, they are even less 
significant. We suppose that the disk material is all contained within —H < z < H, there being 
vaccuum outside this region. Here we and below we find it convenient to take H to be the 
maximum value of the semi-thickness of the disk. 

We begin by finding ^>' Gh . We remark that interior to the disk satisfies the Poisson equation 

V 2 *^ = 4*Gf/ h , (A2) 

while external to the disk it satisfies Laplace's equation. For the warping perturbations we 
consider, p' h and ^' Gh are odd functions of z. Thus d^>' Gh /dz is an even function of z and 
accordingly it has the same values on the top and bottom of the disk. 

To solve ( |A2[) , we use the standard method of Hankel transforms (see, e.g., Binney &: 
Tremaine 1987). That is we write (after removing the factors exp(zy? + at)): 

®Gh = / %(k,z)J rn (kr) kdk, (A3) 
Jo 
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oo 



Ph = / p' k (k,z)J m {kr) kdk, (A4) 



o 



where p' k (k,z) and ^f' k (k,z) are the Hankel transforms of p' h (r,z) and ^' Gh (r,z) respectively. The 
Bessel function is denoted by J m in standard notation. Here m = 1. For the global disturbances 
we are interested in, the significant values of k ~ r _1 at some characteristic radius. 



The Poisson equation ( |A2| ) gives 

- k 2 % = A7rGp' k . (A5) 

Above the disk, where there is no source (p' k = 0), the solution is ^' k (k, z) = C{k) exp(— kz), where 
C{k) is an arbitrary function which depends only on k. Because $f' k is antisymmetric with respect 
to reflection in the midplane, below the disk (z < 0) we have ^' k (k,z) = — C{k) exp(fcz). 

To find C(k), we multiply ( |A5| ) by z and integrate through the thin disk to obtain 

2 [ H ^r ~ *' k ) + ~ k2 J-H z ^' kdz = 4nG J-oo zp ' kdz ' (A6) 

where the subscript + denotes evaluation in the vaccuum at the upper boundary of the disk where 
z = H. 



From (A6) we find to within a multiplicative error of at most of order \kH\, this being a small 



quantity for global perturbations: 

/oo 
zp' k dz. (A7) 
-co 

From this it similarly follows that at z = H, 

_A = 2nGk / zp' k dz. (A8) 

V OZ J _j_ J— oo 

From equation (|A3| ), it similarly follows that at z = H, 



dz 



/oo poo 
/ zp' k J m {kr) k 2 dk dz. (A9) 
-oo J 



We recall that for global density perturbations we expect that p' k is significant only for k ~ r 1 at 
a characteristic radius. 

Returning to the expression of interest given by equation ( |A1| ), we first use Poisson's equation to 
express densities in terms of potentials: 
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G 



4ttG 



(AlO) 



However, we shall retain only the derinatives with respect to z in the Poisson operators which 
leads to a multiplicative error of order H 2 /r 2 for globally varying functions. This in turn leads 
to a multiplicative error of order H/r in the final expression we obtain so it may be neglected. 
Using ( A10| ) in this way in ( |A1| ) we obtain 



pf h dz 



i r H (d 2 ^ G a*' 



Gh 



4ttG J-h \ dz 2 dz 



+ d 2 v Gh d* G \ dz 

dz 2 dz 



1 



2ttG V dz dz 



(All) 



where an integration by parts has been performed and the symmetry properties of the potentials 
have been used. 

We also note that for the equilibrium potential 



G 



dz 



2?rGS. 



(A12) 



Using this together with (A9) in (All), we obtain 



/oo poo poo 

pf' h dz = 2vrGS / / zp' k J m (kr) k 2 dk dz 
-oo J —oo J 



(A13) 



Using ( |A4[ ) this may be equivalently expressed as 



/oo poo 
pf' h dz = 2vrGS / zp' h dz 
-oo J — oo 



S-oo Io° zp' k J m {kr) k 2 dk dz 
I-oo fo° zp' k J m (kr) k dk dz 



(A14) 



We remark that the expression in square brackets gives a mean value of k which for global 
perturbations may be taken ~ r" 1 at some characteristic radius. Further from the expression ( |3l| ) 
we can simply estimate that for global perturbations 



(A15) 



These taken together with the condition Q ~ 1, yields the required estimate 



/oo 
pf' h dz ~ n 2 K ZH%/r 2 . 
-oo 



(A16) 
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Fig. 1.— Toomre parameter Q = kc/(ttGS) for M D (R) = 1.5.10 _3 M and {H/r) max = 1.5.10 -3 . 
The other parameters are given in the text. 
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a) Real part of g and imaginary part of -g for a=0.1 
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b) Real part of g for a=10 5 and a=0.1 
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Fig. 2. — a) Real and imaginary parts of g and —g respectively for a self-gravitating viscous disk 
with a = 0.1. The other parameters are given in the text, b) Real part of g for the same disk but 
a = 10" 5 and a = 0.1. 




Fig. 3. — 3D view of the warped disk, with parameters given in the text. The vertical scale has 
been magnified. The plane z = represents the plane of the undistorted disk. 
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Fig. 4. — The phase angle rj(r), the variation of which defines the degree of twist, plotted in degrees 
as a function of r. 



